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Abstract 


The behavior of an air breathing fuel cell (ABFC) operated on dry-hydrogen in dead-ended mode is studied using theoretical analysis. A 
one-dimensional, non-isothermal, combined heat and mass transport model is developed that captures the coupling between water generation, 
oxygen consumption, self-heating and natural convection at the air breathing cathode. The model is validated against planar ABFC experimental 
measurements over a range of ambient temperatures. The model confirms the strong effect of self-heating on the water balance within passive 
ABFCs. Model analysis provides several conclusions: (1) thermal runaway caused by inadequate heat rejection predominantly limits ABFC 
performance. (2) The natural convection boundary layer represents a significant barrier to cathode mass and heat transfer. (3) Because the mass 
and heat transport numbers associated with natural convection are small, even slight forced convection dramatically affects cell behavior. (4) 
Performance optimization requires maximizing heat rejection while minimizing flooding. Decoupling the latter two phenomena is challenging due 
to the exponential relationship between water vapor saturation and temperature. 
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1. Introduction 


Air breathing polymer electrolyte membrane fuel cells (PEM- 
FCs) use free convection airflow to supply oxygen to their 
cathodes. These cells are typically characterized by low out- 
put power densities compared to forced-convection fuel cells. 
They are nevertheless attractive for portable-power applications 
[1,2] where the simplicity of free-convection oxidant delivery 
can outweigh the cost, complexity, noise and the parasitic power 
consumption introduced by active system design. However, the 
inability to regulate the air stream conditions (flow stochiome- 
try, temperature, humidity) makes water-balanced operation of 
an air breathing fuel cell particularly challenging. 

Planar air breathing cells are characterized by an open cath- 
ode structure that allows ambient air to advect and diffuse to 
the surface, while simultaneously collecting current from a gas 
diffusion layer [3,4]. Planar air breathing fuel cells offer the 
advantage of a flat form factor suitable for small portable appli- 
cations. There have been several recent experimental [3—10] and 
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modeling [11-16] studies on planar air breathing fuel cells. In 
one of the few combined experimental and modeling studies of 
an air breathing PEMFC, Mennola et al. concluded that their 
fuel cell was limited by water removal at 40°C and by oxygen 
depletion at 60°C [14]. In related work, the same group high- 
lighted the importance of developing a significant temperature 
gradient between the fuel cell and the air to drive efficient oxy- 
gen transport to the cathode [17]. Li et al. theoretically analyzed 
convective mass transfer at a free-breathing fuel cell cathode but 
did not consider thermal gradients [13]. Their analysis concluded 
that free-breathing fuel cell performance is significantly limited 
at high current densities due to oxygen mass transfer consider- 
ations. Ying et al. developed a 3D model for an air breathing 
PEMFC that considered both heat and mass transfer [11,12]. 
Their model was applied to channel-type air breathing fuel cells, 
rather than planar open-cathode fuel cells. In such cells, perfor- 
mance was shown to depend strongly on the cathode channel 
geometry. Litster et al. [16] computationally resolved the natu- 
ral convection and mass transfer above the open cathode of an air 
breathing non-planar PEMFC. Their results indicated significant 
self-heating in fuel cells cooled only by natural convection. 

In this paper, we develop a simple one-dimensional engineer- 
ing model for a planar air breathing fuel cell that considers the 
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mixed heat and mass transfer situation at the open cathode. Our 
model is based on a theoretical analysis of free convection heat 
and mass transfer, and these dynamics are coupled to the rate 
of heat production and species generation/depletion at the fuel 
cell cathode. This simplified engineering analysis provides an 
intuitive understanding of the primary factors driving the perfor- 
mance of planar air breathing fuel cells. Furthermore, the model 
qualitatively reproduces the experimentally observed behavior 
of a planar air breathing fuel cell over a range of operating con- 
ditions (Tambient 10-30 °C). Lastly, the model is used to explore 
the limits of air breathing fuel cell operation. 


2. Free convection heat and mass transfer analysis 


The natural convection above the surface of an open fuel cell 
cathode is characterized by the simultaneous transport of ther- 
mal energy, momentum, and chemical species. Heat production 
within the fuel cell creates temperature gradients above the cath- 
ode surface, while oxygen consumption and water production 
give rise to concentration gradients of molecular oxygen and 
water vapor. The air density above the cathode surface, and con- 
sequently the rate of natural convection, is affected by all three 
gradients. In this section, we apply free convection heat and mass 
transfer analysis to determine the heat, oxygen and water transfer 
coefficients (hr, ho, and hy, 0, respectively) that describe trans- 
port within the boundary layer between the cathode surface and 
free air. 

Natural convection driven exclusively by thermal gradients 
is described by the non-dimensional heat transfer coefficient, 
Nuz, the Nusselt number, based on a characteristic system length 
L. Empirical formulae exist (see e.g. ref. [18]) for many prac- 
tical geometries and boundary conditions that correlate Nuz 
with other non-dimensional groups such as the Grashof num- 
ber (Grz) and the Prandtl number (Pr). In general, Nuz = 
f(GrL, Pr). In an analogous fashion, natural convection driven 
exclusively by concentration gradients may be described by the 
non-dimensional mass transfer coefficient, Shy, the Sherwood 
number. The Sherwood number, in turn, can also be related to 
other non-dimensional groups as Shy = f(Grz, Sc), where Sc 
is the Schmidt number. 

Similarity considerations for natural convection driven 
by concurrent temperature and concentration gradients yield 
Nur = f(Grz, Pr, Sc) and Sh; = f(Grz, Sc, Pr). However, 
as a first approximation the simpler correlations given for Nuz 
and Shy may be used to determine the heat and mass transfer 
coefficients provided the effects of temperature and concentra- 
tion on density are properly accounted for [18]. The Prandtl 
number is defined as Pr=v/a, where v is the kinematic vis- 
cosity and « is the thermal diffusivity. In a similar manner, the 
Schmidt number is defined as Sc = v/D, where D is the diffusion 
coefficient. For natural convection driven by concurrent temper- 
ature and concentration gradients the general form of Grashof 
number is: 


_ 8(Ap/p)L* _ 8(s0 = ps)L* 


Gry z z (1) 


Here g is the gravity, L the characteristic dimension of the prob- 
lem, p= (ps + Poo)/2 the average density across the boundary 
layer and A p/p is the density change across the boundary layer 
(which, for a fuel cell cathode, will depend on both the temper- 
ature and concentration changes across the boundary layer). As 
a reasonable approximation, an air breathing fuel cell cathode 
surface can be treated as a heated horizontal plate. The com- 
monly used formula for heat transfer above a heated horizontal 
plate is [18]: 

Nur = 0.54(Grz Pr)'/*, 104 < Gr Pr < 107 (2) 


Based on heat and mass transfer analogy, the mass transfer 


relationship is: 
Shy = 0.54(Gr,Sc)'*, 104 < Gr, Sc < 10’ (3) 


The heat transfer coefficient is therefore: 
-Nurk 2g Ap\ 1 
r= = = o.sau( “4 e) Loa) (4) 


And the oxygen and water vapor mass transfer coefficients are: 


O , 1/4 

Sh 2 Do k: 3/4 Ap 2g = 

Oo = HUL re L a, air = 0.54 D0. air ae a La 0/4 (5) 
a 


HO 1/4 
Sh7?~ Deyo ,ai 3/4 Ap 2g -1/4 
hig = e _ 2" D oa —— | 14°” 


La o 


Note that for the fuel cell cathode model, when calculating the 
heat transfer number, the characteristic length corresponds to 
the entire length of one edge of the fuel cell device (which we 
designate by L). In contrast, the appropriate characteristic length 
involved in H20 and O2 mass transfer is taken as the square root 
of the cathode active area (which we designate by La). 

Accurate fuel cell thermal modeling requires radiation to be 
included in the analysis. At high current densities (when the 
cathode surface reaches 70-80 °C), radiation accounts for 30% 
of the total heat flux from the fuel cell. Failure to account for 
radiation heat transfer leads to large discrepancies between the 
calculated and experimentally measured fuel cell heat balance. 
Modifying Eq. (4) to account for radiation yields: 


hr = hr + hT RAD 
= hy + 2€plateOpottz(T@p, + Tamp)(Topt + Tame) (7) 


where Ar is the modified heat transfer number, AT RAp accounts 
for radiative heat transfer, £plate the plate emissivity of the cath- 
ode and opo}tz is the Stephan—Boltzmann constant. The factor of 
two accounts for radiation heat transfer from both sides of the 
device. 


3. Model development 
3.1. Model scope and assumptions 


The air breathing fuel cell model developed in this paper is 
a simple, one-dimensional, steady state, cathode only fuel cell 
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model. A one-dimensional free-convection model is especially 
well suited to describing the situation at a planar air breath- 
ing fuel cell cathode because there are no directional forced 
flows. The following list details the major assumptions used in 
developing this model: 


One-dimensional transport; 

Steady state conditions; 

Single phase flow (only water vapor transport is considered); 

The effects of liquid water accumulation (flooding) are not 

treated; 

Dead-ended anode; 

Dry (zero humidity) hydrogen supply; 

No water accumulation in the anode; 

No net water transport through the membrane (a ramification 

of flux balance with dead-ended anode assuming no liquid 

water accumulation in anode); 

e Water activity is uniform across the membrane and is in equi- 
librium with the water vapor activity at the cathode catalyst 
layer; 

e The cathode catalyst layer is infinitely thin (treated as a sur- 

face). 


A major limitation of the current model is its failure to account 
for liquid water flooding. Although the model tracks liquid 
water accumulation, the model does not evaluate liquid water 
transport. Also, the model does not take into account the delete- 
rious effects of flooding on fuel cell performance. Experimental 
research suggests that flooding is a significant problem in air 
breathing fuel cells [21]. Therefore, a major emphasis of future 
work will be to implement liquid water transport and flooding 
effects into the model. 

Since water activity is assumed uniform within the membrane 
and no net water transport is permitted to cross it, we assume 
that the conductivity of the membrane is set by the temperature 
and water vapor activity at the cathode catalyst interface. This 
simplification allows us to narrow our model domain so that 
only transport processes occurring between the cathode catalyst 
interface and quiescent ambient air above the cathode surface 
are considered. This region is further broken down into two 
domains: (1) the cathode gas diffusion layer (GDL) and (2) the 
cathode natural convection boundary layer. We apply diffusion 
governing equations to describe the transport processes within 
the GDL and convection governing equations to describe the 
transport processes in the natural convection boundary layer 
just above the cathode surface. The complete model space is 
diagrammed schematically in Fig. 1. Although the equations pre- 
sented in this model development are based on concentration, for 
gas phase species at low pressure, the ideal gas law can be used to 
convert between concentration and pressure as P; = RTc;, where 
P; is the partial pressure of species i, R the ideal gas constant, T 
the temperature and c; is the concentration of species i. 


3.2. Model domain and governing equations 


Space does not permit full discussion of the model devel- 
opment here. Instead, Table 1 summarizes the set of equations 


1D-Model, boundary conditions, governing equations 


catalyst Convection ambient 


—= GDL —— boundary layer —*] 


Ogpr | 
Ż -p =- pim- 
mass |nF ox San. J = hp AC = hy (Come — Cane) | 
heat | jape plo gq =h,AT = hy Lam Toni) | 
Ox Sont. 


Fig. 1. Schematic diagram of the model domain, governing equations and rel- 
evant parameters. Shown are traces of temperature, water concentration and 
oxygen concentration. The relevant flux constant is shown in parentheses by 
each curve. 


used in the model formulation while Table 2 summarizes the 
nomenclature. Model solution is accomplished using an itera- 
tive approach. The boundary values above the cathode natural 
convection boundary layer (Tams and c;,amp) are set by the 
ambient conditions (ambient temperature, pressure and rela- 


Table 1 
Summary of equations used to formulate the air breathing fuel cell combined 
convection and heat transfer cathode model 


Cathode GDL domain governing equations 


de; , oe 
j= -n FD = npp OL ic Where De = = Dj aie (8) 
! dx i GDL ; T 
, dT Topi — Te 
q = —kanL — = —kepi £ (9) 
dx GDL 
Natural convection boundary layer governing equations 
J = hi Aci = hi(Ci AMB — Ci,GDL) (10) 
4 = hr AT = hy(Tamp — Ton) a1) 


Fuel cell heat flux 


T ASRXN,H20@) XH,0,1 A Hcond 
TEGS y ave } cal 12 
q=J] f ocv ) ( 2F ( 2F ) (12) 


Fuel cell polarization model 


V= Eocv — Nkinetic — NIR (13) 
Eocy = FSo3K = 
RT, "05 
| e in | CO2-AMB a4) 
2F CO>.c 


In refs. [24,25] 
nir = İR = jAcen R = jAceu[Relec + Rmem], where Rmem = 


mem 


—— | and oar = 
AcellOa,T 
[3.46a> + 0.016 1a” + 1.45a — 0.175] x [e!268(1/303)-d/D))] (15) 
In ref. [26] 
RT. jco, AMB 
Nkinetic = Ja F In ee ee (16) 
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Table 2 

Summary of model equation symbols and nomenclature 

Symbol and description Units 

Constants 
Universal gas constant, R 8.315Jmol-! K 
Faraday’s constant, F 96,500 C mol`! 
Gravity, g 9.8 ms? 


Stephan-Boltzmann constant, O'poltz 9.67 x 1078 W m~? K4 


Ambient conditions 


Ambient pressure, P 101000 pA 
Ambient relative humidity, RHAMB As specified (%) 
Ambient temperature, TAMB As specified (°C) 


Model terms 


Effective diffusivity, pt mê s7! 


<= 


GDL porosity, ¢ Dimensionless 
GDL tortuosity, T Dimensionless 
Binary diffusivity of species i in air, Dj air m? s7! 
Catalysts layer temperature, Te K 

GDL interface temperature, TgpL K 

GDL thickness, dgpr M 

Thermal conductivity of the GDL, kgpi Wm !K 
Concentration of species i in ambient, cj, AMB mol m~? 
Concentration of species i in GDL, ci.GgpL mol m7? 
Concentration of species i in catalyst layer, cj, mol m~? 
Entropy of fuel cell reaction, ASRXN,H20 (e) —44.43 J mol! K 
Heat of water condensation, AH ona 40kJ mol`! 
Fraction of liquid water, xy,0,) Dimensionless 
Reference reversible fuel cell voltage, Ea K V 

Lumped cell electrical resistance, Relec Q 

Cell membrane resistance, Rmem Q 

Cell active area, Ace m? 

Nafion thickness, dmem m 

Water vapor activity, a Dimensionless 
Charge transfer coefficient, œ Dimensionless 
Reference exchange current density Te K 
Exchange current density activation energy, kJ mol"! 


AGact 


tive humidity). At the catalyst layer, the boundary conditions 
Gj, Te, Cic) are provided as guesses by the solver program for 
a given input voltage (V). Based on these initial guesses, the 
solver algorithm uses Eqs. (8)—-(11) (in Table 1) to re-compute 
ambient conditions. The magnitude and direction of the error 
between the computed ambient conditions (based on the guessed 
inputs for the catalyst layer boundary conditions) and the fixed 
ambient conditions is then used to refine the guess for the cat- 
alyst layer boundary conditions and the solution procedure is 
repeated. Iterative refinement proceeds until the errors between 
the real and calculated ambient conditions drop below a min- 
imum threshold (the error between the input and calculated 
cell voltage is minimized simultaneously as well, using Eqs. 
(13)-(16)). Once acceptable solutions have been calculated at a 
given input voltage, the voltage is stepped to a new value and the 
solution algorithm is repeated. In this manner, full model data 
for a complete IV curve is acquired. The iterative guess refine- 
ment solver algorithm used in this model is based on a form of 
the generalized reduced gradient optimization method that was 
developed by Leon Lasdon, of the University of Texas at Austin, 
and Allan Waren, of Cleveland State University for optimization 
of non-linear problems [19]. 


Table 3 


Example model output from air breathing fuel cell simulation at ambient conditions = 10 °C, 40% RH 


Heat and mass transfer 


GDL surface 


Catalyst surface 


Fuel cell overpotentials 


IV 


Nuz 


RHopi 
(%) 


45 


(atm) 


Pi1,0,GDL 


Po,,GDL 


(atm) 


Topi 


CC) 


(0) 


RH (%) XH,0 


(atm) 


Pr,0,¢ 


(atm) 


nr (V) RQ) T.CO Po,,c 


Nkinetic (V) 


Eocy (V) 


J (mA cm~?) 


Veell (V) 


NA 


6.2 


NA 
6.6 
7.1 


NA 
12.4 


.0 


0.006 
0.014 


0.210 


10.0 


0.00 
0.52 


0.5 


45.0 
100.0 
100.0 


0.006 
0.018 


0.210 


10.0 


90.00 
26.91 


0.00 
0.02 
0.03 
0.03 
0.04 
0.05 
0.06 
0.07 
0.08 
0.09 
0.10 


0.11 


0.00 
0.40 
0.42 
0.44 
0.46 
0.47 
0.49 
0.50 
0.52 
0.53 
0.55 
0.56 
0.59 
0.61 
0.62 
0.63 
0.63 
0.64 
0.64 


1.23 
1.23 
1.23 


6.8 


T2 
Ts. 
8.1 


7.6 
8.1 


4.4 
5.2 
5.9 
6.6 
7.1 


8.5 
8.9 


9.2 


8.4 
8.7 
9.0 


9.3 


9.5 


7.6 
8.1 


9.8 


10.1 9.5 


8.5 
8.9 
9:7 
20.0 
20.1 


9.8 


10.3 


10.2 


10.8 


78.7 


0.200 
0.196 
0.192 
0.188 
0.184 
0.180 


15.8 


0.194 
0.188 
0.182 
0.175 
0.168 
0.160 


0. 


15.9 


76.1 


0.016 


18.4 
21.3 
2 


6 


0.022 


18.6 


26.32 


73.7 


0.019 


71.5 


0.022 


69.5 


0.026 


67.9 


0.032 


66.5 


0.038 


65.3 


046 
0.056 


64.3 


66 


63.5 


0.067 


0.080 62.8 


56 


61.7 


0.115 


10.4 


11.0 


11.1 


51.1 


0.123 
0.127 
0.130 
0.131 


10.5 


43.5 


3 


11.2 10.6 


20.3 


7.7 


11.3 10.7 


20.4 


33.2 
29 


11.3 10.8 


20.4 


T 


0.131 


0.57 
0.5 


00.0 
100.0 


4.6 


6 


28.2 


0.53 


00.0 


0 
0.45 
0. 


00.0 
100.0 


0.171 
0. 


39.6 


40 


43.5 


0.33 
0. 


00.0 
100.0 


0.161 
0. 


47.4 


26 


51.2 
5 


0.19 


0. 


00.0 
100.0 


0.026 
0.031 


0.039 


0.059 


51 


0.074 


0.091 


0.111 


34 


0.135 
0.198 
0.215 


21.6 


25.70 
25.05 


24.39 


:23 
1.23 


25.0 


28.6 


.23 
1:23 


32.4 


23.76 


0.143 
0. 


40.3 


44.4 


0.124 


48.4 


52.4 


0.146 
0.141 
0.139 
0.137 
0.137 
0.136 
0.136 


9.1 


01 


0.093 


60.5 


64.7 


0.00 
0.00 
0.00 
0.00 
0.00 
0.00 


83.0 
70.8 


0.084 


66.3 


69.1 
72.9 


0.223 


0.079 


70.8 


0.076 0.229 61.3 


74.8 
78.2 


23 
1.22 


22.59 


22.06 
21.57 


22. 
1.22 


22. 
1.22 


20.30 


0.12 
0.16 


0. 


24.22 


22 
1.22 
1.21 
1.21 
1.21 


1.21 


28.83 


20 


34.19 


0.24 
0.28 
0.32 
0.37 


76.1 


0.232 54.1 


0.074 
0.073 


40.04 


46.37 


78.9 


48.3 


0.233 


81.1 


1.23 
0.80 
0.78 
0.75 
0.73 
0.70 
0.68 
0.65 
0.63 
0.60 
0.58 
0.55 
0.50 
0.45 
0.40 
0.35 
0.30 
0.25 
0.20 


79 
113 


150 


192 
239 
288 
339 
392 
447 
503 
560 
683 
731 
756 
772 
778 
778 
774 


10.8 


11.4 


26.8 20.5 


0.131 


81.4 


43.6 


0.073 0.233 


83.6 


53.19 


121 


All model parameters are as listed in Table 2. 
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Once a complete solution has been generated, model out- 
puts include calculated IV data, iR and kinetic overpotentials, 
values for the temperature and species concentrations at the cat- 
alyst and GDL surfaces, and the condensed water fraction. In 
analyzing the model results, we find that both water concen- 
tration and relative humidity are useful. To convert between 
RH and water concentration, a standard polynomial expression 
for the saturation pressure of water (Psat,H,0) is used [20], 
yielding: 

Pu 


_“th0 y 100% 
PSAT, H20 


RH = 


R x T x cmo 
1062- 1794+0.02953T—9.1837x 10-9724 1.4454x 1077 T3 


5 X100% 
(17) 


Under certain conditions, the model predicts water vapor 
supersaturation (RH > 100%) at the cathode catalyst layer. Under 
these conditions, liquid water condensation would occur. Our 
model accounts for liquid water condensation with a liquid water 
faction term, x},0,1- 


RH (%) — 100 
100 


The liquid water fraction term represents the fraction of water 
produced at the cathode that is condensed as liquid water. If 
RH < 100%, then xq,0,; = 0. Note that this term appears in the 
fuel cell heat flux Eq. (12) to account for the heat released upon 
liquid water condensation. 

Example output of the model for an air breathing fuel cell 
operating at ambient conditions of 10°C and 40% RH is pro- 
vided in Table 3. As shown in the table, Nuesselt numbers 
obtained from our model are typically in the order of 10-20 
while Sherwood numbers are typically in the order of 6-12. 
These large values indicate the relative dominance of convec- 
tive (rather than diffusive) heat and mass transport above the 
cathode surface. 


XxH,01 = (for RH > 100%) (18) 


4. Model validation and discussion 
4.1. Model validation 


In this section, the model is validated against experi- 
mental data over a range of ambient conditions in order to 
gauge the extent of its applicability. This comparison is made 
against recently described experimental measurements from a 
3cm x 3cm air breathing fuel cell device [21]. The baseline 
properties used by the model are listed in Table 4. Fixed inputs 
and cell properties were obtained from the literature or directly 
from the geometry of the experimental cell. Only the four prop- 
erties highlighted in gray (Kgp i; Relec, ae x and g) were used 
as free-floating fitting parameters. A single set of these fitting 
values was used to fit all measurements and conditions; these 
best-fit values are documented in Table 4. Starting from this 
baseline, the model was validated against measurements made 
at ambient temperatures of 10, 20 and 30°C (in all cases at 


40% RH). A comparison between the experimental and results 
is provided in Fig. 3a—d. 

As Fig. 2 indicates, the model nicely predicts the relative 
trends in air breathing fuel cell behavior as a function of ambient 
temperature. The model recovers the experimentally observed 
decrease in performance with increasing ambient temperature 
due to fuel cell overheating and consequent membrane dry-out. 
This dry-out phenomenon is indicated by the results in Fig. 2c 
and d, which show the effect of cell self-heating (Fig. 2d) and 
the resulting increase in fuel cell resistance above a critical 
current density due to membrane dry-out (Fig. 2c). In partic- 
ular, note how the model effectively captures the rapid increases 
in cell resistance and temperature that occur due to the non- 
linear effects of membrane dry-out at high current densities. 
A traditional non-coupled or isothermal fuel cell model would 
not capture these non-linearities. The model also captures the 
general trends for water accumulation, although not the proper 
magnitude, as shown in Fig. 2b. The model’s superficial treat- 
ment of liquid water likely explains this quantitative discrepancy. 
As discussed in Section 3.1, the model assumes that there is no 
water transport across the membrane, and hence no accumula- 
tion of liquid water in the anode. However, this assumption is 
problematic since anode liquid water accumulation is observed 
experimentally. Experimentally, water condenses in the dead- 
ended anode due to the existence of temperature gradients 
between catalyst layer and anode gas channel walls (which our 
model neglects). Due to this temperature gradient, under condi- 
tions where cathode is close to saturation, hydrogen gas becomes 
supersaturated at the walls of the anode gas channels, leading 
to condensation. Because the anode is dead-ended, any water 
that condenses in the anode will simply accumulate in the anode 
compartment. This additional water sink likely explains why 
experimental measurements of liquid water accumulation are 
greater than the model predictions. 

As shown in Table 3, in addition to predicting overall IV 
behavior, the model provides detailed results for temperature 
and species concentrations at both the catalyst and GDL sur- 
faces. Recently, we have experimentally investigated the actual 
temperature and species concentrations at the GDL surface of 
an air breathing fuel cell by employing microscale oxygen, rel- 
ative humidity and temperature sensors [22]. By comparing the 
experimental and model results for GDL surface temperature and 
species concentration as a function of current density, it is there- 
fore possible to further authenticate the physical correctness of 
our model. An example of this comparison is provided in Fig. 3 
for an air breathing fuel cell operated at ambient conditions of 
21°C, 35% RH. 

Fig. 4 shows that most of the experimentally observed trends 
are also reflected in the model predictions. The temperature fit 
between the experimental and model data is especially good, 
indicating that the heat transfer portion of our model is accurate. 
The model overestimates both Po, and Py,9. The model’s over- 
estimate for Py,0 is likely due to the assumption of zero water 
transport to the anode. As discussed earlier, we experimentally 
observe that some of the water produced by the cell transports to 
the anode. This anode water sink therefore reduces the water con- 
centration at the cathode. The approximate two-fold discrepancy 
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Table 4 
Properties (and their values) used in the base case air breathing fuel cell model 


Property 


Value 


Constants 
Universal gas constant, R 
Faraday’s constant, F 
Gravity, g 
Stephan-Boltzmann constant, Opoiz 
Ambient conditions 
Ambient pressure, P 
Ambient relative humidity, RHamp 
Ambient temperature, T,mp 
Fixed Inputs 
Plate emissivity, plate 
Thermal conductivity of air, Kair 
Binary diffusivity of O2 in air, Do2,sir 
Binary diffusivity of H20 in air, Dino sir 
Reference reversible fuel cell voltage, E’ 303k 
Entropy of fuel cell reaction, Asexn 12012) 
Heat of water condensation, AH ona 
Cell Properties 
Characteristic dimension of fuel cell device (square), L 
Length of active cell side (square), La 
Cell active area, Asen 
Nafion thickness, Smem 
GDL thickness, dep. 
GDL porosity, € 
GDL tortuosity, T 
Exchange current density activation energy, AGaa 
Fitting parameters 
GDL thermal conductivity, kep 
Lumped cell electrical resistance, Rete 
Reference exchange current density, j°303k 


Charge transfer coefficient, a 


8.315 J/mol*K 
96,500 C/mol 

9.8 m/s? 

9.67x10° W/m**K* 


101000 pA 
As specified (%) 
As specified (°C) 


0.90 

0.027 W/m*K 
2.1x10° m/s 
2.6x10° m?/s 
1.23 V 

-44.43 J/mol*K 
40 kJ/mol 


3.0 
50 kJ/mol 


10 W/m*K [27] 
0.012 Q [21] 
5x10° A/cm? 
0.28 [28,29] 


Cell geometry and MEA properties are based on the experimental 3 cm x 3 cm air breathing fuel cell described in ref. [21]. Also see refs. [27-29]. 


between the experimental and model Py,0 readings indicates 
that as much as half of the generated water may transport to the 
anode. This is consistent with the results in Fig. 2b, where the 
model values for accumulated water also deviate by a factor of 
approximately two. 

The largest apparent discrepancy between the experimental 
and model results arises in the humidity data, with the model dra- 
matically overestimating GDL humidity (although adequately 
capturing the shape of the humidity versus current density 
curve). The large discrepancy between the model and measure- 


ments of RH is a direct result of the model’s overprediction of 
Pm,o. Because RH is linearly proportional to Pp,o, the model’s 
overestimate of Py,0 (by approximately a factor of two), leads 
to an overestimate in RH (also by approximately a factor of two). 

Overall, the model capturse the major trends associated with 
cell temperature, dry out and cell-to-atmosphere mass transport, 
while falling short of giving accurate humidity values immedi- 
ately above the fuel cell cathode. We are currently refining our 
model further to account for flooding and water transport to the 
anode. 


124 R. O’Hayre et al. / Journal of Power Sources 167 (2007) 118-129 


—s— Exp: 10C 


1.1 —A— Exp: 20C 
4 —S— Exp: 30C 
. ---m --- Model: 10C 
0.94%. ~~~ Model: 20C 


---@-- - Model: 30C 


voltage [V] 


O 100 200 300 400 500 600 700 800 
current density [mA/cm?] 


~ 
O 
~ 


—B— Exp: 10C 
—A— Exp: 20C 
—O— Exp: 30C 
---W--- Model: 10C 
---Ġ&--- Model: 20C 
---@--- Model: 30C 


resistance [mQ] 


O 100 200 300 400 500 600 700 800 
current density [mA/cm?] 


= 
EE 
=~ 


1.8 


| —8— Exp: 10C 
1.6 —A— Exp: 20C 
1.474 —S— Exp: 30C 
---lll--- Model: 10C 
1.24 ---à--- Model: 20C aN 
---@--- Model: 30C Ea 
14 M- S 
i 7 


Le Le y 
O 100 200 300 400 500 600 700 800 
current density [mA/cm?] 


accumulated water after 2h [g] 


(d) 100 

90 |S Be: 10¢ 
= —Å— Exp: 20C 
o 80 |_4 Exp: 30C 
g 70 4- -W ---Model: 10C 
Z an |å- Model: 20C 
S 60 ---® --- Model: 30C 
a 507 E 
§ 404 

30 


O 100 200 300 400 500 600 700 800 
current density [mA/cm?] 


Fig. 2. Comparison of experimental and model results for an air breathing fuel cell operated at ambient temperatures of 10, 20 and 30°C, 40% relative humidity. (a) 
Polarization curves, (b) accumulated water after 2 h of steady state operation at a range of current densities, (c) cell resistance as a function of current density and (d) 


GDL surface temperature as a function of current density. 
4.2. Scenario modeling 


Having validated the model against experimental air breath- 
ing fuel cell data, we now use the model to investigate several 
scenarios. Specifically, we examine the effect of GDL thickness, 
GDL thermal conductivity and increasing heat/mass transfer 
with the addition of forced convection. In all scenarios, ambient 
conditions of 20°C and 40% RH are chosen to best reflect the 
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Fig. 3. Comparison of experimental (closed symbols) and model (open symbols) 
results for GDL surface relative humidity, temperature and species concentra- 
tions as a function of current density. The model and experimental IV curves are 
also shown for comparison. 


likely conditions encountered by an actual air breathing fuel cell 
device. 


4.2.1. Effect of GDL thickness 

Fig. 4 provides a parametric analysis of how variations in the 
GDL thickness affect modeled air breathing fuel cell behavior. 
Inspection of Fig. 4 reveals that the model predicts an optimum 
in fuel cell performance for a GDL thickness of approximately 
500 um. Thicker or thinner GDLs result in decreased perfor- 
mance. (Please note that this model investigation does not take 
the electrical resistance of the GDL into account; these electrical 
resistance differences may be important when comparing very 
thick versus very thin GDLs.) Interestingly, the real GDL thick- 
ness used in current air breathing fuel cells (~300 um, indicated 
by the base case) is slightly thinner than the predicted model 
optimum. The reason for a model optimum is uncovered from 
careful analysis of the complete model data. For the thinnest 
GDL, water is rejected too easily from the cathode, leading to 
premature membrane dry-out. For the thickest GDL, the thermal 
resistance of the GDL itself starts to become an important fac- 
tor, causing a more rapid temperature rise at the catalyst layer 
and leading to dry-out. This temperature effect is observable 
in the increased slope of the curves in Fig. 4d with increas- 
ing GDL thickness. This scenario indicates that GDL thickness 
optimization involves a delicate balance between the heat and 
water transport rates within the fuel cell. This issue of matching 
the heat rejection and water rejection requirements is one of the 
most important considerations in optimizing air breathing fuel 
cell performance. 
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Fig. 4. The effect of GDL thickness on modeled air breathing fuel cell behavior. The base case (GDL thickness = 300 um) is indicated with solid triangles. (a) 
Polarization curves, (b) accumulated water after 2h of steady state operation at a range of current densities, (c) cell resistance as a function of current density and 


(d) catalyst/membrane temperature as a function of current density. Ambient conditions are 20°C, 40% RH. The other parameters used in this model simulation are 
detailed in Table 4. 
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Fig. 5. The effect of GDL thermal conductivity (kapL ) on modeled air breathing fuel cell behavior. The base case (kap = 10 W mK! is indicated with solid triangles. 
(a) Polarization curves, (b) accumulated water after 2h of steady state operation at a range of current densities, (c) cell resistance as a function of current density and 
(d) catalyst/(membrane temperature as a function of current density. Ambient conditions are 20°C, 40% RH. The other parameters used in this model simulation are 
detailed in Table 4. 
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4.2.2. Effect of GDL thermal conductivity 

Fig. 5 provides a parametric analysis of how variations in 
the GDL thermal conductivity affect modeled air breathing fuel 
cell behavior. With increasing GDL thermal conductivity, air 
breathing fuel cell performance first increases, and then satu- 
rates for thermal conductivity values greater than 30 W mK™!. 
At first, performance is improved by increasing kgp_, which 
improves heat rejection from the cathode and delays the onset 
of membrane dry-out. However, the saturation behavior for 
kapt > 30 W mK~—! indicates that for highly conductive GDLs, 
the thermal bottleneck shifts from the GDL to the convection 
boundary layer. Thus, heat advection away from the GDL by 
natural convection in the free air above the cathode becomes 
the new rate-controlling process. As with the previous scenario, 
it is clear that fuel cell performance is improved by increasing 
the heat rejection capabilities of the cell. By rejecting greater 
amounts of heat, the onset of dry-out can be delayed to higher 
and higher current densities. Our model suggests that it is this 
membrane dry-out process, and not oxygen depletion, that limits 
the performance of air breathing fuel cells. 

In order to operate an air breathing fuel cell at extremely high 
power densities, significantly greater heat rejection capabilities 
than those obtainable from natural convection are required. For 
example, a fuel cell operating at 0.5 V and 2.0 A cm~? (cor- 
responding to a power density of 1 W cm~?) would produce 
approximately 2 W cm7? of heat. To maintain a steady state 
catalyst temperature of 70 °C (approximately the maximum sus- 
tainable temperature before the onset of thermal runaway) would 
require a heat transfer number of approximately 400 W m~? K 
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(assuming Tamp = 20°C). Advection of this heat load from the 
cathode surface requires a Nusselt number of approximately 
1000. Such a Nusselt number is certainly not obtainable with 
natural convection and would correspond to extremely vigor- 
ous forced convection. Rather than relying on convection, one 
interesting alternative would be to remove this heat via conduc- 
tion to a heat sink. Such thermal engineering possibilities will 
be an interesting subject of future experimental and modeling 
activity. 


4.2.3. Effect of forced convection 

To reinforce the conclusions drawn from the previous sce- 
nario, this final scenario examines what type of performance 
benefits may be possible if forced convection is introduced into 
the fuel cell model. As shown in Table 3, the Nusslet/Sherwood 
numbers obtained due to natural convection above the fuel cell 
cathode are in the order of 5—20. This range is on the same order 
of magnitude as naturally occurring air currents in an indoor 
building. In contrast, a light (but perceptible) breeze would yield 
Nusselt/Sherwood numbers of 30-50 while vigorous forced air 
flow can yield Nusselt/Sherwood numbers of 100-1000. Since 
even a light breeze yields higher transport numbers than those 
due to natural convection, it is interesting to investigate what 
effect this might have on air breathing fuel cell performance. This 
analysis is provided in Fig. 6, which summarizes the changes in 
model fuel cell behavior due to the superposition of an addi- 
tional Nusselt/Sherwood term accounting for varying degrees 
of forced convection. Increases of 10, 30, 100 and 300 to the 
Nusselt/Sherewood numbers are considered. As with the previ- 
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Fig. 6. The effect of forced convection on modeled air breathing fuel cell behavior. The base case (natural convection only) is indicated with solid triangles. (a) 
Polarization curves, (b) accumulated water after 2h of steady state operation at a range of current densities, (c) cell resistance as a function of current density and 
(d) catalyst/membrane temperature as a function of current density. Ambient conditions are 20 °C, 40% RH. The other parameters used in this model simulation are 


detailed in Table 4. 
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Fig. 7. Water balance in an air breathing fuel cell. The properties plotted in 
this graph were obtained from model simulation of an air breathing fuel cell 
at 20°C, 40% RH. All other properties are as listed in Table 4. Water genera- 
tion increases linearly with increasing current density. Water removal increases 
exponentially with increasing current density (due to cell self-heating). The 
competition between water generation and removal leads to three character- 
istic operating regions: humidification (at low current density), flooding (at 
intermediate current densities) and dry-out (at high current density). 


ous scenarios, all other parameters remain fixed as specified by 
the base case properties listed in Table 4. 

As indicated by Fig. 6, forced convection dramatically 
improves fuel cell performance by increasing heat transfer 
(through improved advection from the cathode surface), thus 
delaying the onset of membrane dry-out. For the highest forced 
convection rate (Nuz/Shz ~ 300), the cell never reaches the crit- 
ical temperature for dry-out and performance is limited only 
by the iR losses. However, as shown by the data in Fig. 6b, 
liquid water accumulation also increases with increasing forced 
convection, indicating that these cells, in reality, would face mas- 
sive flooding problems and never achieve the predicted model 
performance. The failure of our model to account for flooding 
limits our ability to attempt further quantitative analysis in this 
scenario. However, these results provide an important conclu- 
sion about the coupling between thermal and water management 
in air breathing fuel cells; improved thermal rejection is desir- 
able because it delays the onset of dry-out, but it also leads 
to increased liquid water accumulation (because exponentially 
less water can be removed at lower cell temperatures). Thus, the 
desire to increase heat rejection often conflicts with the desire 
to avoid flooding. This “Goldilocks” problem of air breathing 
fuel cell design will be further discussed in the following section 
(Fig. 7). 


4.2.4. Water balance in air breathing fuel cells 

As indicated above, the difficulty of reconciling heat and 
water management often leads to flooding or dry-out problems in 
air breathing fuel cells. Flooding is most easily tracked by moni- 
toring accumulated water (a signal of flooding) while dry-out can 
be detected by monitoring cell resistance (a dramatic increase in 
cell resistance indicates membrane dry-out). For example, the 
experimental and model results in Fig. 2b and c show the char- 
acteristic hallmarks of both flooding (at intermediate current 
densities) and dry-out (at elevated current densities). Not sur- 


prisingly, operation at higher ambient temperature reduces the 
window of flooding and expedites dry-out. Our model can be 
used to understand the factors that control flooding and dry-out. 
This analysis is presented in Fig. 7, which plots cell resistance 
and water accumulation (our signatures for flooding and dry- 
out) along with the principle water flux curves responsible for 
the water balance in the air breathing fuel cell. 

The changing water balance dynamics in the fuel cell lead 
to the three regimes indicated on the figure: self-humidification 
at low current density, saturation and flooding at moderate cur- 
rent densities and dry-out at high current densities. Ideal water 
balance is only achieved at two operating points: (1) at the tran- 
sition between the self-humidification and flooding regimes and 
(2) at the transition between the flooding and membrane dry 
out regimes. The complex interplay between heat and water 
balance, which leads to these transition points, can be under- 
stood by carefully examining the water flux curves. The regime 
transitions are controlled by the competition between the water 
generation rate and the maximum water removal rate. The gen- 
erated water flux versus current density curve is linear with 
slope = 1/2F. The maximum water flux that can be removed 
from the cathode increases exponentially with increasing cur- 
rent density due to the coupling between current density and 
cathode surface temperature. The cathode surface temperature 
increases roughly linearly with increasing current density (see, 
for example; Fig. 2d), but this leads to an exponentially increas- 
ing water vapor saturation pressure at the cathode surface, and 
therefore leads to an exponential increase in the flux of water 
that can be removed from the cathode. Note that the curves 
intersect twice. At current densities below the first intersec- 
tion, self-humidification occurs (as the two curves converge with 
increasing current density). Between these two intersections, 
the water flux generated by electrochemical reaction exceeds 
the maximum water removal rate, leading to membrane sat- 
uration and eventual flooding. At current densities above the 
second intersection, dry-out occurs (as the two curves diverge 
with increasing current density). 

A change in the ambient conditions (ambient temperature 
or humidity) will cause the water flux removal curve to shift. 
For example, Fig. 8 shows the impact of increasing the ambi- 
ent humidity on water balance. Increasing the ambient humidity 
causes the cathode water removal curve to shift downwards (a 
moisture ambient atmosphere decreases our ability to remove 
water from the cathode.) This shift expands the region of mem- 
brane saturation, intensifying flooding, but delaying the onset of 
dry-out. A change in ambient temperature causes similar effects. 


4.2.5. Properties of the natural convection boundary layer 

Because forced convection can dramatically impact air 
breathing fuel cell behavior, it is insightful to ask: how far does 
the natural convective boundary layer penetrate into the ambi- 
ent air above the fuel cell cathode? The characteristic thickness 
of the natural convection boundary layer can be approximated 
from our model using the following simple relation: 


D 
dcony © 7 (19) 
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Fig. 8. Influence of ambient humidity on air breathing fuel cell water balance. 
Relative to the situation at 40% RH (closed symbols), increasing the ambient 
RH to 90% shifts the cathode water removal curve downwards and exacerbates 
liquid water accumulation (open symbols). Model simulation at Tamp = 20°C, 
all other properties as listed in Table 4. 


The room temperature diffusivities for air and water (see 
Table 4) are Do, =2.1 x 107° m2/s_ and Dyno = 2.6 x 
1075 m?/s, respectively. The model-calculated mass transfer 
numbers for O2 and H20 are approximately ho, = 0.006 m/s 
and hy,9 = 0.007 m/s. Eq. (19) therefore yields an estimated 
thickness of the natural convection boundary layer of approx- 
imately 0.0035-0.0037 m=3.5-3.7 mm. This is an order of 
magnitude larger than the thickness of the typical fuel cell GDL 
(~300 um). This surprisingly large boundary layer thickness 
has recently been confirmed by spatially resolved experimental 
measurements of oxygen, water and temperature profiles above 
an air breathing fuel cell GDL [22]. The large thickness of the 
convective boundary layer is one of the key distinguishing char- 
acteristics of free-breathing fuel cells compared to forced-flow 
fuel cells. In a forced flow situations, the convective boundary 
layer thickness is greatly reduced (to the order of 30-100 um), 
consistent with the one to two orders of magnitude increase 
in h. 

Because the model considers coupled heat and mass trans- 
port, itis instructive to ask what fraction of the natural convection 
cathode mass transfer is driven by the thermal gradient, the water 
vapor gradient and the oxygen gradient, respectively. This sep- 
aration can be accomplished via application of the Boussinesq 
approximation [23], which allows the relative contributions of 
AT, Aco, and Acp,o to Ap/p to be isolated from one another 
(for small density variations). For a fuel cell cathode, assuming 
that density changes are only caused by gradients in tempera- 
ture, water and oxygen, the Boussinesq approximation may be 
written as: 
=e = B" (Topi — Tame) + B°?(co;,GDL — CO;,AMB) 


+ p'?°(cH,0,GDL — CH,0,AMB); 
B™°(cH,0,GDL — CH:0,AMB) < 1, 
B°2(co),GDL — CO2,AMB) 1, B' (Topi — Tame) < 1 
(20) 
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Fig. 9. Relative contribution of the thermal, oxygen, and water vapor gradients 
on the natural convection process at an air breathing fuel cell. Results are based 
on the standard model simulation at ambient conditions of 20°C, 40% RH, all 
other properties as listed in Table 4. 


Here 6T is the thermal volumetric expansion coefficient and 
B°2and 642° are the solutal volumetric expansion coefficients. 
The thermal volumetric expansion coefficient is defined as: 


pa (2) ; (21) 
P 


And the solutal volumetric expansion coefficient is: 


_1 (4 
m, 5 


Applying this analysis to our standard air breathing fuel 
cell model simulation provides the results illustrated in Fig. 9. 
In this figure, the relative contribution from each of the three 
gradients is plotted as a function of current density for a typ- 
ical operating condition (Tamp =20°C, RH=40%). As shown 
by the figure, the thermal gradient represents the largest driv- 
ing force for natural convection, while the water vapor and 
oxygen gradients show considerably less influence. The oxy- 
gen gradient is less important than the water vapor gradient 
because the molar weight difference between Oz and N3 is 
much smaller than molar weight difference between H2O 
and N2. 


5. Conclusions 


The one-dimensional, non-isothermal model presented in this 
paper has been designed to capture the coupling between water 
generation, oxygen consumption, self-heating and natural con- 
vection at the cathode of an air breathing fuel cell. Although 
simple, the model captures the general trends observed in exper- 
imental air breathing fuel cells operated over a wide range of 
ambient temperatures. The model has been used to explore the 
major factors determining air breathing fuel cell performance. 
This analysis has resulted in several important conclusions for 
the fuel cell parameters explored: 


(1) Thermal runaway caused by inadequate heat rejection 
predominantly limits air breathing fuel cell performance. 
Oxygen depletion is a less important factor. 
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(2) In air breathing fuel cells, the natural convection boundary 
layer represents a significant barrier to mass and heat trans- 
fer and is a primary limitation to fuel cell performance. As 
documented in Table 3, over half of the temperature and 
concentration drop occurs in the air layer rather than the 
GDL. Natural convection is therefore a primary limitation 
in air breathing fuel cell performance. 

(3) Because the mass and heat transport numbers associated 
with natural convection are small, even a minor amount of 
forced convection (e.g. due to motions of room air) dramat- 
ically affects fuel cell behavior. 

(4) Optimal air breathing fuel cell design requires maximizing 
heat rejection while also minimizing the effects of flood- 
ing. Decoupling these two phenomena is challenging due to 
the exponential relationship between water vapor saturation 
pressure and temperature. 


Our model analysis of thermal runaway suggests that air 
breathing fuel cells should be designed to operate at membrane 
temperatures below 70 °C in order to avoid catastrophic dry-out. 
Unfortunately, the model also shows that passive air breathing 
fuel cell operation at temperatures below 70 °C is complicated by 
a different water problem: flooding. The contradicting require- 
ments to avoid dry-out but also to avoid flooding reflect the 
“Goldilocks” problem of air breathing fuel cell operation. In 
the ideal case, the rates of water generation and removal should 
match, and not just at one operating point, but across the range of 
operating conditions and current densities. Aligning the water 
generation and removal rates over a wide range of operating 
conditions requires decoupling the water removal rate from its 
exponential temperature dependence. 

Although our model permits tracking of liquid water accumu- 
lation, it currently fails to account for deleterious performance 
losses due to liquid water flooding. As both the model and experi- 
mental results indicate, liquid water accumulation is a significant 
factor in air breathing fuel cell operation at intermediate current 
densities. A major emphasis of future work will be to imple- 
ment anode water transport and flooding effects into the fuel 
cell model. 

The conclusions provided by the combined heat and mass 
transfer air breathing fuel cell cathode model presented in 
this paper have important design ramifications. One intrigu- 
ing possibility is to deploy passive heat/mass transport “fin” 
or “chimney” structures to improve cathode performance by 
magnifying natural convection. Passive designs to mitigate 
flooding (via thermal or gravity assist, for example) can also 
be considered. The fact that even a small degree of forced con- 
vection significantly improves performance suggests that even 
air breathing fuel cells in the 1-5 W power range may benefit 
from an integrated low-power fan or “gas-stirrer” to improve 
convection. 
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